Page 122
% solve the logistic eqn with r=1, K=3
% for starting values y0=1, y0=2, y0=4
% Make up an m-file, fig432.m, as follows
%     function yprime = fig432(t,y)
%     r=1; K=3; yprime=y.*r.*(1-y./K);
tspan=[0 5];
[t1,y1]=ode23('fig432',tspan,1);
[t2,y2]=ode23('fig432',tspan,2);
[t4,y4]=ode23('fig432',tspan,4);
plot(t1,y1,t2,y2,t4,y4)

                   Page 123
K=1; rho=3; c=(rho/K)/(1+rho);
y(1)=0.05; t=[1:60];
for i = 2:60
   y(i)=(1+rho)*y(i-1)*(1-c*y(i-1));
end
plot(t,y)

                   Page 124
r=1; theta=0.2; K=1; y=0:0.05:1;
f=r.*(y/theta - 1).*(1-y/K);
plot(y,f); hold on
xaxis = zeros(size(y));
plot(y,xaxis)

                   Page 125
% Make up an m-file, fig434.m:
%     function yprime = fig434(t,y)
%     % with r=1, theta=.2 , and K=1.
%     r=1; theta=0.2; K=1;
%     yprime = y.*r.*(1-y./K).*(y/theta-1);
tspan=[0 3];
[t05,y05]=ode23('fig434',tspan,.05);
[t1,y1]=ode23('fig434',tspan,.1);
[t3,y3]=ode23('fig434',tspan,.3);
[t5,y5]=ode23('fig434',tspan,.5);
[t7,y7]=ode23('fig434',tspan,.7);
[t15,y15]=ode23('fig434',tspan,1.5);
plot(t05,y05,t1,y1,t3,y3,t5,y5,t7,y7,t15,y15)

                  Exercises/Experiments
1c.
tt=0:10:200;
pop=[3.929214, 5.308483,7.239881,...
   9.638453, 12.866020, 17.069453,...
   23.191876, 31.433321, 39.818449,...
   50.155783, 62.947714, 75.994575,...
   91.972266, 105.710620, 122.775046,...
   131.669275, 151.325798,179.323175,...
   203.302031, 226.545805, 248.709873];
plot(tt,pop,'x'); hold on;
alpha=387.980205; beta=54.0812024;
 delta=0.02270347337;
Anderfit=alpha./(1+beta*exp(-delta*tt));
plot(tt,Anderfit)

3d.
% make an m-file, exer43.m
%   function Bprime=exer43(t,B);
%   r=.48; K=15; a=2; b=2;
%   Bprime=r*B.*(1-B/K)-a*B.^2./(b^2+B.^2);
K=15; a=2; b=2; r=.48;
B=0:.1:20; plot(B,Bprime)
[t,y1]=ode23('exer43',[0 30],1);plot(t,y1)
 hold on
[t,y2]=ode23('exer43',[0 30],2);plot(t,y2)
[t,y4]=ode23('exer43',[0 30],4);plot(t,y4)
[t,y5]=ode23('exer43',[0 30],5);plot(t,y5)
[t,y6]=ode23('exer43',[0 30],6);plot(t,y6)
[t,y8]=ode23('exer43',[0 30],8);plot(t,y8)
[t,y10]=ode23('exer43',[0 30],10);plot(t,y10)
[t,y12]=ode23('exer43',[0 30],12);plot(t,y12)
[t,y14]=ode23('exer43',[0 30],14);plot(t,y14)
[t,y16]=ode23('exer43',[0 30],16);plot(t,y16)

4.
% make an m-file, delayFcn0.m:
%   function y=delayFcn0(t)
%   y = 1;
N=10; % # steps per unit interval
delT=1/N; % so delta t=0.1
% t is now linked to index i by t=-1+(i-1)*delT
% set init. val's via delay fcn f0
for i=1:N+1
   t=-1+(i-1)*delT; f(i)=delayFcn0(t);
end
% work from t=0 in steps of delT
tfinal=10; % end time tfinal, end index is n
% solve tfinal=-1+(n-1)*delT for n
n=(tfinal+1)*N+1;
K=3; r=1.2;
for i=N+1:n-1
   t=-1+(i-1)*delT;
   delY=r*f(i-N)*(1-f(i-N)/K)*delT; % N back=delay of 1
   f(i+1)=f(i)+delY; % Eulers method
end
t=-1:delT:tfinal; plot(t,f);